Chapter 6 Beta diversity

load("data/data.Rdata")
beta_q0n <- genome_counts_filt %>%
  column_to_rownames(., "genome") %>%
  filter(rowSums(. != 0, na.rm = TRUE) > 0) %>%
  dplyr::select_if(~!all(. == 0)) %>%
  hillpair(., q = 0)

beta_q1n <- genome_counts_filt %>%
  column_to_rownames(., "genome") %>%
  filter(rowSums(. != 0, na.rm = TRUE) > 0) %>%
  dplyr::select_if(~!all(. == 0)) %>%
  hillpair(., q = 1)

beta_q1p <- genome_counts_filt %>%
  column_to_rownames(., "genome") %>%
  filter(rowSums(. != 0, na.rm = TRUE) > 0) %>%
  dplyr::select_if(~!all(. == 0)) %>%
  hillpair(., q = 1, tree = genome_tree)

beta_q1f <- genome_counts_filt %>%
  column_to_rownames(., "genome") %>%
  filter(rowSums(. != 0, na.rm = TRUE) > 0) %>%
  dplyr::select_if(~!all(. == 0)) %>%
  hillpair(., q = 1, dist = dist)

6.1 location

6.1.1 Richness diversity

betadisper(beta_q0n$S, sample_metadata$location) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df  Sum Sq  Mean Sq      F N.Perm Pr(>F)    
Groups     5 0.51463 0.102926 11.076    999  0.001 ***
Residuals 86 0.79917 0.009293                         
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
               Aruba     Brazil  CaboVerde      Spain    Denmark Malaysia
Aruba                1.0000e-03 3.0000e-03 1.0000e-03 1.0000e-03    0.012
Brazil    8.1559e-05            9.6000e-02 7.7000e-02 9.1000e-01    0.132
CaboVerde 6.2922e-04 1.1254e-01            1.0000e-03 2.7000e-02    0.941
Spain     1.7363e-12 7.3505e-02 3.5940e-05            2.7000e-02    0.002
Denmark   1.3206e-07 9.0556e-01 2.9720e-02 2.3526e-02               0.049
Malaysia  9.6128e-03 1.3730e-01 9.2726e-01 4.6570e-04 5.9273e-02         
adonis2(beta_q0n$S ~ location, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q0n$S))), 
        permutations = 999) %>%
        broom::tidy() %>%
        tt()
tinytable_haqvhsrb2vfhkkqgn2im
term df SumOfSqs R2 statistic p.value
location 5 6.57607 0.2442751 5.559605 0.001
Residual 86 20.34468 0.7557249 NA NA
Total 91 26.92075 1.0000000 NA NA
pairwise.adonis(beta_q0n$S, sample_metadata$location, perm = 999)
                   pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1        Aruba vs Brazil  1 2.0520609 7.280171 0.19528265   0.001      0.015   .
2     Aruba vs CaboVerde  1 1.0026264 3.213494 0.09975614   0.001      0.015   .
3       Aruba vs Denmark  1 1.9513181 7.018855 0.19486614   0.001      0.015   .
4      Aruba vs Malaysia  1 1.1650328 3.678697 0.11257172   0.001      0.015   .
5         Aruba vs Spain  1 2.0433431 8.167535 0.21974917   0.001      0.015   .
6    Brazil vs CaboVerde  1 2.1090004 9.164611 0.24013374   0.001      0.015   .
7      Brazil vs Denmark  1 0.5703370 2.907951 0.09113563   0.002      0.030   .
8     Brazil vs Malaysia  1 0.9283278 3.953419 0.11996993   0.001      0.015   .
9        Brazil vs Spain  1 0.6343957 3.769479 0.11503016   0.001      0.015   .
10  CaboVerde vs Denmark  1 1.8099512 8.070070 0.22373314   0.001      0.015   .
11 CaboVerde vs Malaysia  1 1.2143821 4.593887 0.14094322   0.001      0.015   .
12    CaboVerde vs Spain  1 1.8141430 9.281722 0.24896172   0.001      0.015   .
13   Denmark vs Malaysia  1 1.0124763 4.418610 0.13629855   0.001      0.015   .
14      Denmark vs Spain  1 0.5305346 3.310771 0.10573905   0.001      0.015   .
15     Malaysia vs Spain  1 0.8450695 4.218747 0.13094076   0.001      0.015   .
#pdf("figures/beta_q0_loca.pdf",width=9, height=5)
beta_q0n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = "sample") %>%
  group_by(location) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = location, fill = location)) +
    geom_point(size = 4) +
    #   stat_ellipse(aes(color = beta_q1n_nmds$Groups))+
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    scale_color_manual(values = location_colors)+
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

#dev.off()

6.1.2 Neutral diversity

betadisper(beta_q1n$S, sample_metadata$location) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq  Mean Sq      F N.Perm Pr(>F)   
Groups     5 0.2728 0.054560 3.7079    999  0.005 **
Residuals 86 1.2654 0.014715                        
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
               Aruba     Brazil  CaboVerde      Spain    Denmark Malaysia
Aruba                9.0000e-03 3.0000e-03 1.0000e-03 1.0000e-03    0.001
Brazil    9.2383e-03            7.4800e-01 2.0500e-01 6.6700e-01    0.878
CaboVerde 7.7521e-04 7.7033e-01            3.0600e-01 8.4800e-01    0.902
Spain     1.6646e-05 2.2327e-01 3.0424e-01            3.9000e-01    0.274
Denmark   3.3254e-04 6.6322e-01 8.7642e-01 3.7291e-01               0.773
Malaysia  1.1422e-03 8.5989e-01 8.9763e-01 2.4754e-01 7.7435e-01         
adonis2(beta_q1n$S ~ location, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))), 
        permutations = 999) %>%
        broom::tidy() %>%
        tt()
tinytable_o3o0iba2cvfctnzpidje
term df SumOfSqs R2 statistic p.value
location 5 5.709562 0.221304 4.88821 0.001
Residual 86 20.090068 0.778696 NA NA
Total 91 25.799630 1.000000 NA NA
pairwise.adonis(beta_q1n$S, sample_metadata$location, perm = 999)
                   pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1        Aruba vs Brazil  1 1.3739907 4.789154 0.13766227   0.001      0.015   .
2     Aruba vs CaboVerde  1 1.0636589 3.806728 0.11603497   0.001      0.015   .
3       Aruba vs Denmark  1 1.6037782 5.813780 0.16699651   0.001      0.015   .
4      Aruba vs Malaysia  1 1.2584980 4.466640 0.13346545   0.001      0.015   .
5         Aruba vs Spain  1 1.4382789 5.558030 0.16083181   0.001      0.015   .
6    Brazil vs CaboVerde  1 1.5901673 7.027154 0.19505161   0.001      0.015   .
7      Brazil vs Denmark  1 0.6954159 3.122217 0.09719804   0.003      0.045   .
8     Brazil vs Malaysia  1 0.7313843 3.199010 0.09935119   0.001      0.015   .
9        Brazil vs Spain  1 0.4007935 1.948927 0.06297236   0.031      0.465    
10  CaboVerde vs Denmark  1 1.8198218 8.556140 0.23405480   0.001      0.015   .
11 CaboVerde vs Malaysia  1 1.3171418 6.019861 0.17695136   0.001      0.015   .
12    CaboVerde vs Spain  1 1.5415454 7.905417 0.22017337   0.001      0.015   .
13   Denmark vs Malaysia  1 1.0691578 4.970172 0.15074753   0.001      0.015   .
14      Denmark vs Spain  1 0.5766622 3.014217 0.09718823   0.006      0.090    
15     Malaysia vs Spain  1 0.6446972 3.265586 0.10444666   0.001      0.015   .
#pdf("figures/beta_q1n_loca.pdf",width=9, height=5)
beta_q1n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = "sample") %>%
  group_by(location) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = location, fill = location)) +
    geom_point(size = 4) +
    #   stat_ellipse(aes(color = beta_q1n_nmds$Groups))+
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    scale_color_manual(values = location_colors)+
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

#dev.off()

6.1.3 Phylogenetic diversity

betadisper(beta_q1p$S, sample_metadata$location) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq Mean Sq      F N.Perm Pr(>F)  
Groups     5 0.2152 0.04304 2.1047    999  0.073 .
Residuals 86 1.7587 0.02045                       
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
               Aruba     Brazil  CaboVerde      Spain    Denmark Malaysia
Aruba                0.06400000 0.05500000 0.00100000 0.09900000    0.002
Brazil    0.06236518            0.87200000 0.33000000 0.66400000    0.669
CaboVerde 0.05279004 0.86912737            0.20400000 0.76400000    0.532
Spain     0.00060645 0.32391452 0.20693940            0.10200000    0.517
Denmark   0.09811567 0.63811338 0.73976492 0.10108967               0.309
Malaysia  0.00519311 0.67845753 0.52055144 0.50696000 0.31212601         
adonis2(beta_q1p$S ~ location, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1p$S))), 
        permutations = 999) %>%
        broom::tidy() %>%
        tt()
tinytable_nwvt1dmxuhi2q76nk369
term df SumOfSqs R2 statistic p.value
location 5 2.690791 0.2683828 6.309563 0.001
Residual 86 7.335152 0.7316172 NA NA
Total 91 10.025944 1.0000000 NA NA
pairwise.adonis(beta_q1p$S, sample_metadata$location, perm = 999)
                   pairs Df  SumsOfSqs    F.Model         R2 p.value p.adjusted sig
1        Aruba vs Brazil  1 0.64714035  6.1133065 0.16928127   0.001      0.015   .
2     Aruba vs CaboVerde  1 0.56229778  5.3504131 0.15575979   0.001      0.015   .
3       Aruba vs Denmark  1 1.03758657  9.5591181 0.24790811   0.001      0.015   .
4      Aruba vs Malaysia  1 0.48047063  5.0352414 0.14794199   0.002      0.030   .
5         Aruba vs Spain  1 0.60047832  6.7919744 0.18976250   0.001      0.015   .
6    Brazil vs CaboVerde  1 0.70961443  8.1023431 0.21837821   0.001      0.015   .
7      Brazil vs Denmark  1 0.32742702  3.5968626 0.11034383   0.022      0.330    
8     Brazil vs Malaysia  1 0.19570506  2.5119802 0.07971509   0.026      0.390    
9        Brazil vs Spain  1 0.06456161  0.9106384 0.03044530   0.459      1.000    
10  CaboVerde vs Denmark  1 1.07120011 11.9405505 0.29895808   0.001      0.015   .
11 CaboVerde vs Malaysia  1 0.50032095  6.5728056 0.19011490   0.002      0.030   .
12    CaboVerde vs Spain  1 0.82104026 11.9236999 0.29866220   0.001      0.015   .
13   Denmark vs Malaysia  1 0.54258755  6.8084757 0.19559821   0.001      0.015   .
14      Denmark vs Spain  1 0.24626102  3.3999417 0.10827860   0.033      0.495    
15     Malaysia vs Spain  1 0.26932060  4.5771932 0.14050299   0.001      0.015   .
#pdf("figures/beta_q1p_loca.pdf",width=9, height=5)
beta_q1n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = "sample") %>%
  group_by(location) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = location, fill = location)) +
    geom_point(size = 4) +
    #   stat_ellipse(aes(color = beta_q1n_nmds$Groups))+
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    scale_color_manual(values = location_colors)+
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

#dev.off()

6.1.4 Functional diversity

betadisper(beta_q1f$S, sample_metadata$location) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq  Mean Sq      F N.Perm Pr(>F)    
Groups     5 1.5016 0.300325 9.1732    999  0.001 ***
Residuals 86 2.8156 0.032739                         
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
               Aruba     Brazil  CaboVerde      Spain    Denmark Malaysia
Aruba                2.0000e-03 6.0000e-01 2.0000e-03 1.0000e-03    0.494
Brazil    8.6955e-04            1.0000e-03 6.7300e-01 2.5400e-01    0.003
CaboVerde 5.8936e-01 9.7727e-04            3.0000e-03 1.0000e-03    0.879
Spain     1.5125e-03 6.8358e-01 1.5804e-03            4.6000e-02    0.004
Denmark   1.7539e-04 2.3957e-01 9.2505e-05 3.3658e-02               0.001
Malaysia  4.9615e-01 1.0071e-03 8.8292e-01 1.5959e-03 7.4748e-05         
adonis2(beta_q1f$S ~ location, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1f$S))), 
        permutations = 999) %>%
        broom::tidy() %>%
        tt()
tinytable_ap1s72pyl0acpd333eg6
term df SumOfSqs R2 statistic p.value
location 5 8.815599 0.5844961 24.19552 0.001
Residual 86 6.266793 0.4155039 NA NA
Total 91 15.082391 1.0000000 NA NA
pairwise.adonis(beta_q1f$S, sample_metadata$location, perm = 999)
                   pairs Df   SumsOfSqs    F.Model           R2 p.value p.adjusted sig
1        Aruba vs Brazil  1  3.77109355 40.6261637  0.575228238   0.001      0.015   .
2     Aruba vs CaboVerde  1  0.81599223  5.5801053  0.161367505   0.023      0.345    
3       Aruba vs Denmark  1  4.54372380 50.8948430  0.637022880   0.001      0.015   .
4      Aruba vs Malaysia  1  0.87672195  6.0397240  0.172367910   0.010      0.150    
5         Aruba vs Spain  1  3.49955224 36.9984243  0.560595570   0.001      0.015   .
6    Brazil vs CaboVerde  1  1.79573937 27.7696018  0.489163230   0.001      0.015   .
7      Brazil vs Denmark  1  0.27026868 35.0545409  0.547260825   0.001      0.015   .
8     Brazil vs Malaysia  1  1.69391038 26.6369446  0.478763612   0.001      0.015   .
9        Brazil vs Spain  1 -0.01258568 -0.9666596 -0.034482498   0.943      1.000    
10  CaboVerde vs Denmark  1  2.85656348 47.6207158  0.629731090   0.001      0.015   .
11 CaboVerde vs Malaysia  1 -0.02164970 -0.1836835 -0.006603446   0.986      1.000    
12    CaboVerde vs Spain  1  1.56212224 23.8545897  0.460028512   0.001      0.015   .
13   Denmark vs Malaysia  1  2.82610262 48.0023359  0.631590271   0.001      0.015   .
14      Denmark vs Spain  1  0.38029019 58.5480612  0.676480332   0.001      0.015   .
15     Malaysia vs Spain  1  1.47580731 22.9256154  0.450178466   0.001      0.015   .
#pdf("figures/beta_q1f_loca.pdf",width=9, height=5)
beta_q1f$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = "sample") %>%
  group_by(location) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = location, fill = location)) +
    geom_point(size = 4) +
    #   stat_ellipse(aes(color = beta_q1n_nmds$Groups))+
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    scale_color_manual(values = location_colors)+
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

#dev.off()

6.2 Behaviour

6.2.1 Richness diversity

betadisper(beta_q0n$S, sample_metadata$origin) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df  Sum Sq   Mean Sq      F N.Perm Pr(>F)
Groups     1 0.00411 0.0041147 0.2442    999  0.626
Residuals 90 1.51634 0.0168483                     

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
         Domestic Feral
Domestic          0.634
Feral     0.62238      
adonis2(beta_q0n$S ~ origin, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q0n$S))), 
        permutations = 999,
        strata = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))) %>% pull(location)) %>%
        broom::tidy() %>%
        tt()
tinytable_jwqzi3yc0w8puuhjp1j7
term df SumOfSqs R2 statistic p.value
origin 1 0.4407951 0.0163738 1.498173 0.178
Residual 90 26.4799587 0.9836262 NA NA
Total 91 26.9207538 1.0000000 NA NA
pairwise.adonis(beta_q0n$S, sample_metadata$origin, perm = 999)
              pairs Df SumsOfSqs  F.Model        R2 p.value p.adjusted sig
1 Domestic vs Feral  1 0.4407951 1.498173 0.0163738   0.078      0.078    
beta_q0n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = join_by(sample == sample)) %>%
  group_by(origin) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = origin, fill = origin)) +
    geom_point(size = 4) +
    scale_color_manual(values = origin_colors) +
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

6.2.2 Neutral diversity

betadisper(beta_q1n$S, sample_metadata$origin) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq   Mean Sq      F N.Perm Pr(>F)
Groups     1 0.0080 0.0080048 0.5581    999  0.432
Residuals 90 1.2909 0.0143435                     

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
         Domestic Feral
Domestic          0.435
Feral     0.45698      
adonis2(beta_q1n$S ~ origin, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))), 
        permutations = 999,
        strata = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))) %>% pull(location)) %>%
        broom::tidy() %>%
        tt()
tinytable_uu0ifuqco0081afeh6eh
term df SumOfSqs R2 statistic p.value
origin 1 0.386668 0.01498735 1.369385 0.314
Residual 90 25.412962 0.98501265 NA NA
Total 91 25.799630 1.00000000 NA NA
pairwise.adonis(beta_q1n$S, sample_metadata$origin, perm = 999)
              pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1 Domestic vs Feral  1  0.386668 1.369385 0.01498735   0.158      0.158    
pdf("figures/beta_q1n_behaviour.pdf",width=9, height=5)
beta_q1n$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = join_by(sample == sample)) %>%
  group_by(origin) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = origin, fill = origin)) +
    geom_point(size = 4) +
    scale_color_manual(values = origin_colors) +
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 
dev.off()
quartz_off_screen 
                2 

6.2.3 Phylogenetic diversity

betadisper(beta_q1p$S, sample_metadata$origin) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df  Sum Sq  Mean Sq      F N.Perm Pr(>F)
Groups     1 0.01229 0.012295 0.5442    999  0.469
Residuals 90 2.03337 0.022593                     

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
         Domestic Feral
Domestic          0.462
Feral     0.46262      
adonis2(beta_q1p$S ~ origin, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1p$S))), 
        permutations = 999,
        strata = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))) %>% pull(location)) %>%
        broom::tidy() %>%
        tt()
tinytable_w7r655ibvm1y4p1ydgdo
term df SumOfSqs R2 statistic p.value
origin 1 0.2191073 0.02185403 2.010807 0.209
Residual 90 9.8068363 0.97814597 NA NA
Total 91 10.0259435 1.00000000 NA NA
pairwise.adonis(beta_q1p$S, sample_metadata$origin, perm = 999)
              pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1 Domestic vs Feral  1 0.2191073 2.010807 0.02185403   0.077      0.077    
beta_q1p$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = join_by(sample == sample)) %>%
  group_by(origin) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = origin, fill = origin)) +
    geom_point(size = 4) +
    scale_color_manual(values = origin_colors) +
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    ) 

6.2.4 Functional diversity

betadisper(beta_q1f$S, sample_metadata$origin) %>% permutest(., pairwise=TRUE) 

Permutation test for homogeneity of multivariate dispersions
Permutation: free
Number of permutations: 999

Response: Distances
          Df Sum Sq  Mean Sq      F N.Perm Pr(>F)  
Groups     1 0.2486 0.248604 3.3502    999  0.055 .
Residuals 90 6.6786 0.074206                       
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Pairwise comparisons:
(Observed p-value below diagonal, permuted p-value above diagonal)
         Domestic Feral
Domestic          0.062
Feral    0.070508      
adonis2(beta_q1f$S ~ origin, 
        data = sample_metadata %>% arrange(match(sample,labels(beta_q1f$S))), 
        permutations = 999,
        strata = sample_metadata %>% arrange(match(sample,labels(beta_q1n$S))) %>% pull(location)) %>%
        broom::tidy() %>%
        tt()
tinytable_ae61qk00t5d0t7on8n65
term df SumOfSqs R2 statistic p.value
origin 1 0.3816603 0.02530502 2.336579 0.185
Residual 90 14.7007312 0.97469498 NA NA
Total 91 15.0823914 1.00000000 NA NA
pairwise.adonis(beta_q1f$S, sample_metadata$origin, perm = 999)
              pairs Df SumsOfSqs  F.Model         R2 p.value p.adjusted sig
1 Domestic vs Feral  1 0.3816603 2.336579 0.02530502   0.113      0.113    
beta_q1f$S %>%
  vegan::metaMDS(., trymax = 500, k = 2, trace=0) %>%
  vegan::scores() %>%
  as_tibble(., rownames = "sample") %>%
  left_join(sample_metadata, by = join_by(sample == sample)) %>%
  group_by(origin) %>%
  mutate(x_cen = mean(NMDS1, na.rm = TRUE)) %>%
  mutate(y_cen = mean(NMDS2, na.rm = TRUE)) %>%
  ungroup() %>%
  ggplot(aes(x = NMDS1, y = NMDS2, color = origin, fill = origin)) +
    geom_point(size = 4) +
    scale_color_manual(values = origin_colors) +
    geom_segment(aes(x = x_cen, y = y_cen, xend = NMDS1, yend = NMDS2), alpha = 0.9, show.legend = FALSE) +
    theme_classic() +
    theme(
      axis.text.x = element_text(size = 12),
      axis.text.y = element_text(size = 12),
      axis.title = element_text(size = 20, face = "bold"),
      axis.text = element_text(face = "bold", size = 18),
      panel.background = element_blank(),
      axis.line = element_line(size = 0.5, linetype = "solid", colour = "black"),
      legend.text = element_text(size = 16),
      legend.title = element_text(size = 18),
      legend.position = "right", legend.box = "vertical"
    )